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SUMMARY 

This paper is concerned with both the solution sensitivity and solution 
accuracy of large dynamic problems involving NASTRAN SOLUTION 8 (ioO., the 
steady state dynamic response option wherein all response quantities vary as 
eiojt^ where to is the driving frequency and t is time). Using a submerged 
steel plate with a viscoelastic layer as the bench mark sample, the solu- 
tion sensitivity and solution accuracy is checked. The solution sensitivity 
is examined by running the same finite element model on different computers, 
different versions of NASTRAN, and different precision levels. The solution 
accuracy is evaluated for these same runs by comparing the NASTRAN results 
with the exact solution of the same problem. 

SYMBOLS 

[B] Damping Matrix 

Cl Dilational Wave Speed in Fluid 

{F} Applied Force Vector 

[K] Stiffness Matrix 

[1<] Modified Complex Stiffness Matrix 

[M] Mass Matrix 

k Wave Number (u/ci) 

P. Incident Fluid Pressure 

Pq Plane Wave Amplitude 

Pg Back Side Fluid Pressure 

Pg Front Side Fluid Pressure (Scattered Component) 

t Time 

{U} Solution Displacement, Vector 
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SYMBOLS (Conf d) 


X 





r 

»p 

i 

>y 


p 


Spatial Coordinate 
Driving Frequency 
Residual Solution Vector 
Real Elastic Lame' Constants 
Corresponding Viscoelastic Constants 
Material Mass Density 


INTRODUCTION 


This paper is concerned with the solution accuracy of 1, 2, or 3-dimen- 
sional steady state (time harmonic) structural and/or continuum problems whose 
response quantities all vary in time in proportion to jhe linear equa- 

tions of motion for such problems usually reduce to an expression of the form 

[-a)2[M ] + im[Bj + [ K]]{U> = {F> , (1) 

[K] 

where [M], [B], and [K] denote the mass, damping and stiffness matrices (MDD, 
BDD and KDD using usual NASTRAN DMAP notation), u is the driving frequency and 
{F} are the applied forces. The results presented in this paper focus on con- 
tinuum type (e.g., figure 1) applications with structural damping, however, 
once the form of equation (1) has been constructed, the solution becomes a 
matter of solving large banded symmetrical systems of complex linear algebraic 
simultaneous equations. Clearly, such equations can also be the end point 
resulting from many other NASTRAN steady state formulations, either from direct 
structural formulations or from related fields through analogies. Thus compari 
sons of solution accuracy, run time, etc. can be viewed and interpreted in a 
more general vein than simply applying only to problems of the type depicted in 
figure 1. 

The motivation for this comparative study resulted as a consequence of 
obtaining some unexpected results on some solution 8 (steady state time har- 
monic rigid format) problems similar to the one shown in figure 2, except for 
the fact that the initial model had inclusions throughout the rubber thus mak- 
ing analytical solutions to the problem unwieldy. 



PARAMETRIC STUDY MODEL 


In order to better understand the accuracy limitations of the results of 
the initially more complicated inclusion filled model, a simpler homogeneous 
layered model (figure 2) was constructed and physically corresponds to a 
totally submerged 2.0" steel plate with a 3.05" viscoelastic rubber layer 
glued to the steel surface. The input corresponds to an incident pressure 
wave 

Pi = Pq k = o)/Ci (2) 

where x is the horizontal coordinate along the line of propagation, c^ is 
the dilatational wave speed in the fluid, and Pq is the plane wave amplitude. 
The exact analytical solution to this problem is known (ref, 1), consequently 
an accuracy check on the finite element solutions is available. Clearly, the 
figure 1 model is a spatially one-dimensional problem, consequently the cor- 
responding finite element model need only be one element wide as was done, 
for example in ref, 1, However, the finite model was made up to eight ele- 
ments wide for the following reasons: (1) the model simulates the more com- 

plex model except for the fact that the inclusions are removed by filling 
their space are with uniform elements having the same material property as 
the surrounding rubber material; (2) the problem is artificially made mathe- 
matically larger so that more meaningful comparisons of CPU run times could 
be made; (3) larger problem sizes tend to draw out any potential problems 
with equation solvers. It is not our intent to discuss or explain the setup 
of wave propagation problems of the type represented by the figures 1-2 
example model; the reader is referred to refs, 1 and 2 for supplementary 
details. In fact, the demonstration problem used here is very similar to 
the one used in the ref, 1 sample problem except that the plate and visco- 
elastic thicknesses are different, the damping coefficient in the visco- 
elastic layer is different and that the steel plate is represented here 
approximately with CBAR elements rather than with solid elements as in ref, 1 
Specifically, the material constants employed are listed below 

DEMONSTRATION PROBLEM PHYSICAL CONSTANTS 


MATERIAL 
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Water 

345;600, 

0.0 

0,0 

0.0 

.000096 

Steel 

17,307,000. 

11,538,000. 

.0,0 

,0.0 

.000735 

Viscoelastic 

Material 

86,703, 

115,9 

8670,3 

11.59 


.0003599 


where the meaning of the elastic and viscoelastic constants are defined in 
detail in ref, 1, 
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Since the topic of interest here is related to the class of problem 
treated by ref» 1, it appears appropriate to print an errata to the ref. 1 
paper 

• in equation (15) of ref. replace G22 = with G22 = x'^ + 2^^^ 

• in equation (16) of ref. 1, replace G22 = e with G22 = e(l + 2pVx^) 

• in equation (2) of ref. 1, replace with 

• in equation (17) of ref. 1, replace +nj with +irij in the l <2 defi- 
nition 2 ^ 


PARAMETER VARIATIONS 


The basic finite element model, figure 2, was exercised for a frequency 
sweep of 7 different incident frequencies (3.0 kHz, 4.0 kHz, 6.0 kHz, 8.0 
kHz, 17.5 kHz, 22,5 kHz, 35.0 kHz). Running the figure 1 model on NASTRAN 
for the above frequency sweep is designated as a typical run and correspond- 
ingly assign it a "run number", which runs from the number 2 through 9, Run 
number 1 is the exact solution and therefore is the only non NASTRAN designa- 
tion.(it is called a run since even the analytical solution involves a computer evalua- 
tion). 

Next, the same frequency sweep input data was rerun while varying the 
following parameters: 

• solution precision (S.P. or D.P. on the same computer) 

• type of computer (UNIVAC 1108; DEC-VAX; CDC Cyber 175) 

• level of NASTRAN (both NASA and MSC versions are considered) 

• date (i.e., the same input is resubmitted on the same computer, 
using the same version of NASTRAN but on different days) 

The last parameter (i.e., the date) seems a waste of computer time, however 
as is shown later, some unexpected results are encountered. 


DMAP INSTRUCTIONS FOR PRINTING SOLUTION ERROR RESIDUAL 


It is of interest to know the accuracy of the solution solving capa- 
bility of the equation solver used by the particular version of NASTRAN 
employed by the user. Specifically, if the solution {U} is found by NASTRAN, 
how well does it satisfy the linear simultaneous equations (1)? Consider 
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substituting the solution {U} into equation (1) and then transposing the 
applied force vector to the left hand side of equation (1) to obtain 

[K]{U> - {F} = {A} (3) 

If the equations have been solved exactly, then the residual vector {a} will 
be identically zero. The appearance of large nonzero entries in {a} would 
imply potential inaccuracies in the solution vector {U}. The question of 
"how large is large?" should be viewed by comparing the size of a particular 
entry in the {A} vector to the size of the applied loads (for this reason, 
the load vector {F} is also printed). For example, a residual of .2 would 
be a big residual if the applied forces are on the order of 1,0 lbs,: however, 
if applied forces are say 100,000 lbs,, the ,2 residual is acceptable. 


In order to print out the residual vector {A} for the 3,000 Hz driving 
frequency case, the following DMAP instructions were used (note, the {A} 
vector is printed with the heading DELS0L), 

• For UNI VAC 1108, SOL 8, LEVEL 17,0 (Runs 3a, 4a) 


ALTER 159 

ADD5 KDD, BDD, MOD,, /KBARX/C,Y,ALPHA=(1,0,0,0)/C,Y,BETA=(0, 0,18849. 5592)/ 

C,Y,GAMA=(-355305758,44,0.0) $ 

MPYAD KBARX,UDVF,PDF/DELS0L/C,N,O/C,N,1/C,N,-1/ $ 

MATPRN DELS0L,PDF,,,// $ 

ENDALTER 

CEND 

where frequency dependent input constants BETA and GAMA are to be i nput by 
the user and are simply defined as: 

BETA = 0.0 + ioj 
GAMA = -0)2 + i 0.0 

where w = the driving frequency in radians/sec 

• For UNIVAC 1108, SOL 8, NASA LEVEL 15.5 (Run 2) replace ALTER 159 
with ALTER 139 

• For VAX, SOL 8, NASA LEVEL 17,5 (Run 5b) same as 1108, NASA LEVEL 
17,0 

• For CDC CYBER, SOL 8, MSC LEVEL 48B (Run 6) replace ALTER 159 with 
ALTER 139 

• For VAX, SOL 8, MSC LEVEL 52 (Run 7) replace ALTER 159 with ALTER 
139 

• For VAX, SOL 8, MSC LEVEL 60 (Run 8) replace ALTER 159 with ALTER 
139 

• For VAX, SOL 26, MSC LEVEL 60 (Runs 9b, 9c, 9d) 
replace ALTER 159 with ALTER 409 

replace UDVF with UHV (3rd line) 
replace PDF with PD (3rd and 4th line) 
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When a large residual is encountered, it is desirable to know the node 
and component number where a large residual appears (i«e», knowing the run 
number of the questionable residual, what node-component number does this 
correspond to?)„ By inserting a DIAG 22 card, the desired correspondence 
between run number of {A} and the node-component number can be made.. 

It is important to note that the simple DMAP sequence as presented will 
apply to only a single frequency; thus, if a frequency sweep is employed, 
only the nw column of the DELSOL vector (i,e., {A}) will be correct; the 
remaining columns of DELSOL should be ignored, where n = the valoe in the 
frequency sweep appearing on the FREQ card» -ard« 


DISCUSSION OF RESULTS 


The primary variables of interest to us in this study are: (1) the 

transmitted pressure in the fluid on the back side of the steel plate, Pg, 
(e,g., in element number 100352 as shown in figure 2) and (2) the scattered 
pressure in the fluid on the front side of the plate, Ps (e.g», in element 
number 100378 as shown in figure 2 ), The transmitted pressure is read 
directly from the NASTRAN printout, whereas the scattered pressure is obtained 
indirectly from the NASTRAN printout by simply subtracting the incident 
pressure (equation (2)) from the total pressure printed by NASTRAN. The 
scattered pressure is of prime importance with regard to establishing the 
energy absorbing properties of the viscoelastic configuration. As discussed 
in ref, 2, it has been our experience that for steady state wave propagation 
problems of the type considered here, at least 10 elements per wave length 
are needed to adequately compute the pressure response for elements of the 
type employed in this study. In order to demonstrate this accuracy limita- 
tion, the model has been purposely exercised in a driving frequency range 
that is too high for the mesh to properly produce sufficiently accurate 
results (i,e,, the mesh is too coarse for some of the higher frequencies). 

For a rubber wave speed of C 2 = 15540, in/sec, and the coarsest element in 
the rubber mesh (,1" x ,1" elements), the following frequency vs, element/ 
(wave length) chart is constructed: 
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freq. (kHz) 

elements/wave length 
(rubber in figure 2) 

3.0 

51.8 

4.0 

38.8 

6,0 

25.9 

8.0 

19.4 

17.5 

8.8 

22.5 

6.9 

35.0 

4.4 


Based on the above chart, it is expected that the accuracy of the finite 
element solution in relation to the exact solution should start to drift at 
frequencies of 17=5 kHz and higher. 

A total of 14 computer runs were made (designated as runs number 2, 3a, 
3b, o o o o 9d) and are tabulated in Table 1 (for the transmitted back side 
pressure Pg) and in Table 2 (for the scattered front side pressure, Ps)„ 

The pressure results are normalized by the magnitude of the incident pressure 
wave, Pqo In addition to the pressure result. Table 1 has two additional 
pieces of information, namely the magnitude of the largest complex residual, 
IaI, appearing in the residual vector {a} as computed by equation (3), (the 
|a| value is obtained by scanning the DMAP printed DELS0L printout and seek- 
ing out the largest absolute value of all the rows of the {a} vector corres 
sponding to the 3,000 Hz frequency case. Note that the residual for only 
the 3,000 Hz case is reported. Also listed is the total CPU time required 
to execute the full frequency sweep solution for the run in question. 


Word Length (Precision) Sensitivity 


The earlier Level 15.5 version of NASA NASTRAN for a UNIVAC 1108 com- 
puter, does not efficiently solve complex systems of equations of the type 
given by equation (1), (i.e., steady state time harmonic rigid format 8) 
when the double precision option is used. Experience on large problems 
(e.g., the size of the one in ref. 2, pg 435) has demonstrated that in the 
solving of the equation (1), NASTRAN has spent literally hours in the decom- 
position operation. In order to obtain reasonable run times, a PARAM 
DEC0M0PT4 card is added to the bulk data in order to force the decomposition 
to work in single precision. In comparing the NASTRAN Run 2 (Table 1) 
(1108-S.P.) results to the exact solution over the frequency range 
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(3. - 8. kHz) where the mesh is sufficiently fine,* it is noted that at 6.0 
kHz, a 154.0% error in the transmitted pressure is experienced. The corre- 
sponding error on the scattered pressure (Run 2, Table 2) is not as severe, 
namely 16.5%. It is noted that the 154% error is the situation that moti- 
vated this entire comparative study. The percent errors at 3. kHz, 4.0 kHz 
and 8 kHz are also much larger than should be expected for the mesh size 
employed. The largest residual, 1 a|, in the residual vector {A}, is .272 at 
3.0 kHz; this is in comparison to a load vector component of the size (.20) 

(the largest residual was not at a loaded node however). This is a further 
indication that the solution resulting from decomposition in single precision 
is not accurate enough. Upon running the same problem on a Level 17.0 version 
of NASA-NASTRAN on an 1108 computer (single precision must be invoked with 
DMAP) in single precision. Run 3a still results in a similar bad solution with 
a similar worse residual [a]. The fact that Level 17.0 uses a different 
decomposition algorithm did not improve the bad results. However, again run- 
ning the same problem on Level 17,0 on the UNI VAC 1108 in double precision 
(the default situation) resulted in excellent results in comparison to the 
exact solution. For example, at the 6.0 kHz, the percent error reduced from 
154.0% down to 0.4% error. Similarly good results were obtained in the entire 
(3.0 - 8.0 kHz range) for both the transmitted and scattered pressure. The 
double precision gave a very small worse residual at 3,000 Hz, |a| = 2.099x10“®, 
which suggests that equation (1) has been solved accurately. In comparing CPU 
times between Run 3a and Run 4a in Table 1, it appears there is little penalty 
in CPU time between single and double precision runs for Level 17.0. Conse- 
quently, there is really no incentive (from a time saving point of view) to 
make Sol. 8 type runs in single precision, as there was for Level 15.5 NASTRAN. 

Again running the same problem on Level 48B, MSC version of NASTRAN on a 
CDC CYBER 175 computer (Run 6), very good results were obtained in relation 
to the exact solution over the (3,0 - 8.0 kHz range); further, excellent con- 
sistency with the UNIVAC 1108 double precision runs is demonstrated at all fre- 
quencies by comparing Run 4a with Run 6. The worst residual, |a[, on the CDC 
computer is not as good as the 1108 double precision run, but this is expected 
since the CDC single precision word length is slightly smaller than the 1108 
double precision word length; however, it should be noted that differences in 
the decomposition algorithms could also account for differences in the worse 
residual, even if the word lengths were the same. 

The Level 52, MSC version of NASTRAN on a DEC-VAX computer (Run 7) gave 
comparible results to the Level 17.0 NASA NASTRAN double precision 1108, and to 
Level 48B, MSC NASTRAN on a CDC-CYBER 175. The word length in double precision 
on the VAX is slightly more than the single precision CDC and slightly less 
than the double precision 1108. It is noted the Run 7 gave the smallest worse 
residual. 


* It should be emphasized that in comparing any NASTRAN result to the exact 
solution, for the purpose of measuring the formulation quality, it should be. 
done only in the frequency range of 3.0 - 8.0 kHz where the mesh satisfies the 
10 element/wave length criterion. 
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Computer Type/Level of NASTRAN Sensitivity 


It is inconvenient to run the same exact version of NASTRAN on different 
computers due to leasing restri cti ons ,, consequently thi s combination was not 
done. Thus running NASTRAN on different computers always involved running a 
different NASTRAN version as well. In scanning the results of Table 1 and 
Table 2, all the runs performed with decomposition precision using word lengths 
between 60 - 72 bits (i.e., Runs 3a, 4a i 6, 7) gave both accurate results in 
comparison to the exact solution (3.0 - 8.0 kHz range) and consistent results 
from machine-to-machine and version-to-version. We have purposely not commented 
on the accuracy of Runs 5, 8 and 9 involving the DEC-VAX computer due to some 
reservations we have regarding the operating conditions of the particular VAX 
on which these runs were made, and will be discussed next. 


Solution Repeatability 


During the process of preparing this collection of comparative runs, an 
unexplained phenomenom (which is still unexplained as of this writing) occurred, 
namely the fact that Level 60 MSC-NASTRAN run on a DEC-VAX computer gave dif- 
ferent results to the same problem upon rerunning the same data. The Run 9 
series of runs were made on rigid format 26 which is comparable (there are dif- 
ferences in the decomposition routine) to NASA-NASTRAN rigid format 8, For 
example the input producing Run 9a was resubmitted over again (producing Run 
9b) so that the residual vector {A} could be printed (employing the DMAP 
instructions given earlier). In comparing the solutions, the results were 
slightly different (e.g., 5,2% at 3.0 kHz). The same input data was again 
rerun on successive days producing Runs 9c and 9d. Run 9c is the closest to 
the more stable results made on the 1108 and CDC computers. 

Using the same DEC-VAX computer facility, the base case input was resub- 
mitted again employing the Level 17.5 NASA/ GODDARD NASTRAN and thus producing 
Runs 5a and 5b. Again the nonrepeatability of the solution on the VAX was 
experienced, this time with an entirely different version of NASTRAN. 

An MSC Level 60 DEC-VAX computer run was made similar to Run 9a through 
Run 9d, except that the older rigid format 8 instead of the newer MSC rigid 
format 26 was employed. The results were poor in comparison to the exact solu- 
tion, further, the worse residual of [a] = .2127 was unacceptably high. No 
reruns of the same input were made on this version and level of NASTRAN. 

In order to demonstrate that repeatable results are possible (a notion we 
usually assume is true on most modern computers), the base case data was rerun 
on the UNIVAC 1108 computer; single precision Runs 3a and 3b were totally 
repeatable as well as double precision Runs 4a and 4b. 

The facts that (1) the VAX computer resulted in nonrepeatable results 
employing two separate versions of NASTRAN and (2) the Level 60 MSC NASTRAN for 
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the VAX computer, Sol. 8 (Run 8) gave poor results, strongly suggests a prob- 
lem with the particular VAX computer on which the runs were made. The follow- 
ing list provides some possible reasons for nonrepeatability: 

• computer central processing drops a bit in the main memory or 
operating register 

• main memory itself drops a bit between storing and retrieving data 

• floating point accelerator drops a bit during calculations 

• disk subsystem (drive or interface drops a bit) 

• computer temperature rises due to air conditioning not keeping 
up with thermal load during the summer months when the runs were 
made; this could result in dropping a bit by one of the four 
above mentioned possibilities. 

Finally it is noted that as of this writing, non NASTRAN users of the same VAX 
computer that made the runs reported here, did not report any repeatability 
problems with computer program results totally unrelated to NASTRAN. 


CONCLUSIONS 


The paper is concerned with both the accuracy of equation solvers and 
with the accuracy of the problem formulation for large dynamic steady state 
problems (e.g., rigid format 8). Based on a series of computer runs on dif- 
ferent versions of NASTRAN on different computers, the following set of con- 
clusions are drawn: 

• NASTRAN solution decomposition algorithms employing less than a 60 
bit word could lead to serious errors in the results (e.g., employ- 
ing 36 bit single precision words on an 1108 computer gave up to 
154% error in the solution with both Level 15.5 and Level 17.0 ver- 
sions of NASA NASTRAN. 

• Correlation between results run on the 1108 double precision Level 
17.0 NASA NASTRAN; CDC-CYBER 175 single precision Level 48B MSC- 
NASTRAN; and DEC-VAX double precision Level 52 MSC-NASTRAN were 
excellent. 

• Unexplained unrepeatability of results were experienced on the DEC- 
VAX computer for both MSC-Level 60 NASTRAN and NASA-GODDARD Level 
17.5 NASTRAN; Level 17.0 of NASA NASTRAN had no repeatability prob- 
lems on an 1108 computer. 

• A minimum of 10 elements per wave length should be used to model 
traveling wave propagation problems of the type treated in this 
paper; coarser meshes lead to increasingly bad results when comparing 
NASTRAN results to the exact continuum solution to the same problem. 
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As a final comment regarding repeatability, it is noted that it is not being 
suggested that this problem is one to be found in all VAX computers employing 
NASTRAN. The spirit of the NASTRAN Colloquium is to share USER'S experiences, 
thus it was felt that our problem should be brought to the attention of the 
NASTRAN user's community in the event , that similar problems are encountered by 
others. 
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FIGURE 1 - TEST PROBLEM 



FIGURE 2 - FINITE ELEMENT MODEL 



• 856 Elements (CQDMEM, (TRMEM, CBAR) 

• 588 Grid Points 


element 100378 





SoPe (Single Precision) 
D,Po (Double Precision) 


TABLE 1 - NASTRAN COMPARATIVE SOLUTIONS 
(BACK SIDE PRESSURE, element 100352) 




SoPo (Single Precision) 
DoP, (Double Precision) 


TABLE 2 - NASTRAN COMPARATIVE SOLUTIONS 
(FRONT SIDE PRESSURE, element 100378) 


' FREQUENCY 

FRONT SIDE SCATTERED PRESSURE, Ps/Po (amplitude) 


DIIM 

SOLUTION SOURCE 

35. kHz 

22.5 kHz 

17.5 kHz 

8, kHz 

6. kHz 

4. kHz 

3. kHz 

NUMBER 

EXACT ANALYTICAL (Ref. 1) 

.04253 

.08580 

.13310 

.34059 

.44266 

.59075 

.59907 

1 

Level 15.5 NASA-NASTRAN, UN I VAC 1008 
Computer, 36 Bit S.P. Decomposition (Sol .8) 

.15093 

.05609 

.11772 

.32372 

.36949 

.59695 

.62678 

2 

Level 17.0 NASA-NASTRAN, UNIVAC 1108 

.15093 

.05603 

.11759 

.32875 

.39233 

,58218 

.60932 

3a 

Computer, 36 Bit S.P. Decomposition (Sol. 8) 

.15093 

.05603 

.11759 

.32875 

.39233 

.58218 

.60932 

3b 

Level 17.0 NASA-NASTRAN, UNIVAC 1108 

.15093 

.05592 

.11744 

.33419 

.44176 

.59134 

.59943 

4a 

Computer, 72 Bit D.P. Decomposition (Sol .8) 

.15093 

; .05592 

.11744 

, .33419 

.44176 

.59134 

.59943 

4b 

Level 17.5 NASA/GODDARD NASTRAN, DEC-VAX 

.14389 

.06016 

.12378 

.29474 

.40715 

.58416 

.59813 

5a 

Computer, 64 Bit D.P, Decomposition (Sol ,8) 

.14385 

.06717 

.11805 

.24716 

.32194 

.50502 

.56985 

5b 

Level 48B MSC-NASTRAN CDC CYBER 175 
Computer, 60 Bit S.P. Decomposition (Sol. 8) 

.15093 

.05591 

.11744 1 

.33419 

.44176 

.59134 

,59942 

6 

Level 52 MSC-NASTRAN, DEC-VAX Computer, 
64 Bit D.P. Decomposition (Sol. 8) 

.15093 

.05592 

.11744 

.33419 

.44176 

.59134 

.59942 

7 

Level 60 MSC-NASTRAN, DEC-VAX Computer, 
64 Bit D.P. Decomposition (Sol. 8) 

.16173 

.05591 

.11745 

.33203 

.43973 

.59726 

.59972 

8 

Level 60 MSC-NASTRAN, DEC-VAX Computer, 

.15304 

.05422 

.11353 

.28611 

.44209 

.55800 

.50161 

9a 

64 Bit D.P. Decomposition (Sol. 26) 

.15100 

.05838 

.11587 

.33329 

.44178 

.59291 

.60261 

9b 

(4 runs on same data) 

.15095 

.05592 

.11822 

.33563 

.44131 

,59211 

.59944 

9c 


.15093 

.05592 

.12031 

.33385 

.45400 

.59269 

.59942 

9d 



